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ABSTRACT 


Fully  nonlinear  mesoscale  model  simulations  are  used  to  investigate  the  momentum  fluxes  of 
gravity  waves  that  emerge  at  a  “far- held”  height  of  6  km  from  steady  unsheared  flow  over  both 
an  axisymmetric  and  elliptical  obstacle  for  nondimensional  mountain  heights  hm  =  Fr~^  in 
the  range  0.1-5,  where  Fr  is  the  surface  Fronde  number.  Fourier-  and  Hilbert-transform  di¬ 
agnostics  of  model  output  yield  local  estimates  of  phase-averaged  momentum  flux,  while  area 
integrals  of  momentum  flux  quantify  the  amount  of  surface  pressure  drag  that  translates  into 
far-fleld  gravity  waves,  referred  to  here  as  the  “wave  drag”  component.  Estimates  of  surface 
and  wave  drag  are  compared  to  parameterization  predictions  and  theory.  Surface  dynamics 
transition  from  linear  to  high-drag  (wave-breaking)  states  at  critical  inverse  Fronde  numbers 
Fr~^  predicted  to  within  10%  by  the  transform  relations  of  Smith  (1989b).  Wave  drag  peaks 
at  Fr~^  <  hm  ^  2,  where  for  the  elliptical  obstacle  both  surface  and  wave  drag  vacillate 
due  to  cyclical  buildup  and  breakdown  of  waves.  For  the  axisymmetric  obstacle,  this  occurs 
only  at  hm  =  1.2.  At  hm  ^2-3  vacillation  abates  and  normalized  pressure  drag  assumes  a 
common  normalized  form  for  both  obstacles  that  varies  approximately  as  Wave  drag 

in  this  range  asymptotes  to  a  constant  absolute  value  that,  despite  its  theoretical  shortcom¬ 
ings,  is  predicted  to  within  10-40%  by  an  analytical  relation  based  on  linear  clipped-obstacle 
drag  for  a  Sheppard-based  prediction  of  dividing  streamline  height.  Constant  wave  drag  at 
hm  ~2-5  arises  despite  large  variations  with  hm  in  the  three-dimensional  morphology  of  the 
local  wave  momentum  fluxes.  Specific  implications  of  these  results  to  the  parameterization 
of  subgrid-scale  orographic  drag  in  global  climate  and  weather  models  are  discussed. 
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1.  Introduction 


Finite  computing  resources  force  global  numerical  weather  and  climate  prediction  models 
to  operate  at  spatial  resolutions  that  prevent  an  explicit  simulation  of  the  full  spectrum  of 
dynamics  that  can  arise  when  surface  flow  impinges  on  mesoscale  orography.  Such  dynamics 
exert  drag  forces  on  the  larger  scale  flow  and  so  all  credible  numerical  weather  and  climate 
prediction  models  incorporate  parameterizations  of  these  unresolved,  or  subgrid- scale,  oro¬ 
graphic  mesoscale  drag  forces  (hereafter  OMD).  Guidance  for  OMD  parameterization  comes 
from  observations,  numerical  modeling  and  laboratory  experiments. 

Mesoscale  dynamics  observed  in  and  around  realistic  three-dimensional  mountain  ranges 
remain  challenging  to  model  and  reproduce  even  with  state-of-the-art  tools  (e.g..  Smith 
et  ah  2006).  More  reproducible  OMD  hndings  have  emerged  from  idealized  numerical  and 
laboratory  simulations  in  which  a  prescribed  upstream  flow  of  constant  x-directed  velocity 
U  and  buoyancy  frequency  N  impinges  upon  an  isolated  three-dimensional  elliptical  obstacle 
h{x,y)  of  peak  height  hm-  A  typical  form,  also  used  in  this  work,  is 

h{x,y)  =  +  +  (1) 

To  hrst  order,  the  dynamical  response  depends  on  the  normalized  peak  obstacle  height 

=  Nhm/U  =  Fr-\  (2) 

where  Fr  is  surface  Fronde  number.  When  hm  1  flow  passes  across  the  obstacle,  gener¬ 
ating  a  linear  orographic  gravity  wave  response  whose  three-dimensional  form  and  areally- 
integrated  momentum  flux  Dl  can  be  computed  analytically  (Phillips  1984;  Smith  1988). 
As  hm  — >  1,  total  pressure  drag  across  the  mountain,  Dp,  begins  to  exceed  the  linear  solution 
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Dl  by  factors  of  up  to  ~2.  These  so-called  “high-drag”  states  form  as  flow  stagnates  aloft 
(wave  breaking)  and  on  the  upstream  slope  (flow  splitting),  triggering  leeside  collapse  of 
isentropes,  fast  downslope  winds,  turbulence  and  downstream  vortices  (Smolarkiewicz  and 
Rotunno  1989;  Miranda  and  James  1992;  Olafsson  and  Bougeault  1996;  Schar  and  Durran 
1997).  At  hm  >2-3  wave  breaking  eases  or  abates  and  pressure  drag  eventually  drops  below 
Dl.  As  fim  increases,  a  progressively  deeper  surface  flow  layer  stagnates  upstream  and  splits, 
so  that  only  the  flow  above  this  layer  can  pass  across  the  obstacle  peak,  generating  gravity 
waves  of  reduced  amplitude  and  wavelength  compared  to  linear  predictions  (e.g..  Ding  et  al. 
2003;  Lindeman  et  al.  2008). 

The  responses  also  vary  with  the  obstacle’s  aspect  ratio 


(3  =  b/a. 


(3) 


Variations  with  (3  in  the  linear  [hm  1)  hydrostatic  limit  were  derived  analytically  by 
Phillips  (1984).  Dynamical  responses  close  to  the  mountain  as  a  function  of  both  hm  and 
f3  were  hrst  proposed  by  Smith  (1989b)  based  on  limits  of  the  linear  hydrostatic  transform 
relations  of  Smith  (1988)  for  onset  of  wave  breaking  [hm  >  hcb)  and  flow  splitting  [hm  >  h^s)- 
Figure  1  reproduces  the  resulting  regime  diagram,  which  demarks  linear  regions  at  hm  < 
he  in  white  from  various  gray-shaded  nonlinear  responses  at  hm  >  he,  based  on  a  critical 
normalized  obstacle  height  (or  critical  inverse  Fronde  number)  he  =  Fr~^  =  min  hcb,hf 
At  /J  3>  1  it  predicts  transition  from  linear  wave  solutions  at  hm  -C  1  to  a  wavebreaking-only 
regime  at  hm  values  slightly  less  than  unity,  prior  to  the  onset  of  both  wave  breaking  and 
flow  splitting  at  hm  >  1.  As  ellipticity  is  reduced  [f3  — >■  I"*")  the  wavebreaking  region  narrows 
and  eventually  disappears  at  ~  1,  due  to  progressively  greater  wave-held  dispersion.  For 
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/9  <  1  a  transition  from  linear  waves  to  pure  flow-splitting  dynamics  occurs  at  progressively 
larger  he  values  as  (3  decreases,  along  with  progressively  smaller  pressure  drags  (Phillips 
1984).  These  predicted  regimes  have  been  tested  and  broadly  conhrmed  using  numerical 
models  and  tank  experiments  (e.g.,  Miranda  and  James  1992;  Olafsson  and  Bougeault  1996; 
Vosper  et  al.  1999;  Bauer  et  ah  2000;  Ding  et  ah  2003). 

Since  Figure  1  was  formulated  by  Smith  (1989b)  for  responses  in  just  the  lowest  few 
kilometers  above  the  mountain,  numerical  and  laboratory  studies  have  focused  on  quantifying 
the  total  pressure  drag  across  the  obstacle  and  its  variation  with  hm  and  (3  and  have  not 
quantihed  the  fraction  of  the  pressure  drag  due  to  gravity  waves  that  propagate  away  from 
the  obstacle  at  larger  (nonlinear)  values.  OMD  parameterizations  require  estimates  of 
both  the  total  pressure  drag  near  the  mountain  and  the  far- held  gravity  wave  component. 
Some  studies  have  investigated  wave  properties  and/or  momentum  huxes  for  specihe  cases 
(Miranda  and  James  1992;  Vosper  et  al.  1999;  Vosper  2000;  Ding  et  al.  2003;  Lindeman 
et  al.  2008;  Wells  et  al.  2008).  Yet,  despite  extensive  literature  on  linear  wave  properties  at 
hm  -C  1  and  high-drag  wavebreaking  states  at  ~  1,  we  are  not  aware  of  any  work  that  has 
as  yet  systematically  documented  both  surface  pressure  drag  and  gravity  wave  momentum 
huxes  from  idealized  three-dimensional  obstacles  over  the  full  range  of  potential  hm  values: 
closest  perhaps  is  the  work  of  Welch  et  al.  (2001)  for  two-dimensional  sinusoidal  topography 
(see  also  Aguilar  and  Sutherland  2006). 

The  resulting  uncertainty  directly  impacts  the  current  generation  of  OMD  parameteriza¬ 
tions,  which  have  been  heavily  inhuenced  by  this  literature  on  drag  due  to  idealized  three- 
dimensional  obstacles.  Most  modern  OMD  schemes  ht  subgridscale  mesoscale  orography 
within  grid  boxes  to  equivalent  elliptical  obstacle  shapes  in  order  to  utilize  their  documented 
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drag  properties  (Lott  and  Miller  1997;  Scinocca  and  McFarlane  2000),  Yet  these  schemes 
must  quantify  both  the  total  surface  drag  and  the  propagating  gravity- wave  component, 
since  gravity  waves  propagate  vertically  and  can  deposit  momentum  flux  at  higher  altitudes, 
while  the  remaining  drag  is  confined  and  applied  near  the  surface.  Given  the  lack  of  quanti¬ 
tative  guidance  on  orographic  gravity  wave  momentum  fluxes  from  the  modeling  literature 
to  date,  simple  theoretical  estimation  methods  are  employed  that  vary  among  OMD  im¬ 
plementations  (Lott  and  Miller  1997;  Gregory  et  ah  1998;  Scinocca  and  McFarlane  2000; 
Webster  et  ah  2003)  and  require  independent  research  results  for  objective  validation. 

Here  we  conduct  fully  nonlinear  mesoscale  model  simulations  of  flow  across  obstacles  of 
the  form  (1)  at  (3  values  of  1  and  3  for  a  range  of  values  shown  with  solid  circles  in 
Figure  1.  Our  focus  is  on  quantifying  the  momentum  fluxes  of  propagating  gravity  waves 
that  emerge  at  some  altitude  Zi  above  the  nonlinear  surface  dynamics  near  the  obstacle,  and 
how  these  “far  field”  wave  fluxes  relate  both  to  the  total  surface  pressure  drag  and  to  simpler 
predictions  used  in  OMD  parameterizations. 


2.  Modeling  and  Diagnostic  Tools 

a.  WRF  Simulations 

We  use  version  3.0  of  the  Advanced  Research  WRF  (Weather  Research  and  Forecast¬ 
ing)  model  (Skamarock  et  al.  2008).  The  WRF  model  is  nonlinear  and  fully  compressible, 
discretized  with  hfth/third-order  hnite  differences  for  the  horizontal/vertical  advection,  and 
third-order  Runga-Kutta  for  the  time  stepping. 


5 


Our  three-dimensional  computational  domain,  {x,y,z),  is  square  horizontally  with  the 
obstacle  h{x,y)  in  (1)  located  slightly  upstream  of  the  domain  center  to  allow  more  room 
for  dynamics  downstream  of  the  obstacle  to  evolve,  such  that  — L/3  <  x  <  2L/3  and 
— L/2  <  y  <  -h/2,  L  =  900  km  and  z  =  0-30  km.  Inflow/outflow  conditions  are  applied 
at  the  lateral  boundaries.  The  horizontal  gridpoint  resolution  [Ax,  Ay)  is  1.5x1. 5  km^. 
Our  0-30  km  vertical  domain  is  divided  into  200  levels  using  WRF’s  pressure-based  vertical 
coordinate,  120  of  which  are  below  2;  =15  km  where  the  vertical  resolution  Az  125  m, 
with  Az  gradually  increasing  above  15  km. 

In  every  experiment  we  adopt  constant  vertical  prohles  of  iV  =  0.01  rad  s“^  and  x-directed 
horizontal  velocity  U  =  -|-10  m  s“^,  then  vary  and  (3  from  experiment  to  experiment.  We 
choose  a  =  10  km  to  minimize  the  rotational  and  nonhydrostatic  linear  drag  modihcations 
that  become  signihcant  at  larger  and  smaller  widths,  respectively  (Klemp  and  Lilly  1980;  Mi¬ 
randa  and  James  1992),  then  set  h  =  (3a.  We  adopt  a  frictionless  (free-slip)  lower  boundary, 
a  1.5-order  turbulent  closure  scheme,  and  omit  Coriolis  effects  (/  =  0)  and  moisture. 

To  minimize  artihcial  reflection  of  gravity  waves  from  the  model’s  rigid  lid,  the  top  half 
of  the  domain  {z  =  15-30  km)  contains  implicit  vertically-varying  Rayleigh  friction  of  the 
Klemp  et  ah  (2008)  form.  Details  are  given  in  section  2d. 

h.  Diagnostic  Fourier  Transform  (FT)  Method 

To  help  characterize  three-dimensional  gravity  wave  helds  in  WRF  simulations,  we  em¬ 
ploy  a  diagnostic  FT  method  introduced  by  Lindeman  et  ah  (2008).  The  method  utilizes 
the  steady  state  (t  — >■  cxo)  form  of  the  time-dependent  Fourier-ray  algorithm  of  Broutman 
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et  al.  (2006).  For  the  case  of  vertical  velocity 


w{x,y,z)  = 


w{k,l,z)  , 


(4) 


—  OO 


the  method  solves  for  w{k,  I,  z)  using  ray  methods  (Broutman  et  al.  2002).  For  the  specihc 


case  of  constant  U  and  N  prohles  considered  here,  the  Fourier-ray  solutions  at  heights  2;  and 


Zi  are  related  as 


(5) 


where  vertical  wavenumber  m  follows  from  the  dispersion  relation 


m  =  — 


(6) 


(fc,  1)  is  horizontal  wavenumber,  and,  since  we  assume  stationary  waves,  intrinsic  frequency 
uj  =  —kU.  To  evaluate  (4)  we  set  the  Fourier  coefficients  with  fc  >  0  equal  to  zero  (assuming 


f/  >  0),  multiply  the  Fourier  coefficients  with  fc  <  0  by  two,  and  take  the  real  part  of  the 


resulting  Fourier  transform  w{x,y,z),  i.e.  w  =  Re(w).  The  choices  of  fc  <  0  and  m  <  0 
ensure  upward  (downward)  group  (phase)  propagation.  Additional  terms  in  the  dispersion 
relation  can  be  included  (Broutman  et  al.  2009)  but  are  omitted  here  since  they  have  minor 
effects.  We  also  omit  trapping  and  tunneling  terms  (Broutman  et  al.  2003,  2009)  due  to  the 
absence  of  shear. 

For  linear  problems  {hm  1)  the  FT  method  can  be  initialized  at  the  surface  Zi  =  0 
using  the  orographic  elevations  h{x,y).  The  lower  boundary  condition 


w{k,  I,  Zi  =  0)  =  —zixh{k,  /), 


(7) 


where  h{k,  1)  is  the  Fourier  transform  of  h{x,  y),  yields  the  linear  steady-state  mountain  wave 
response  using  (4)- (6). 
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In  characterizing  mountain  waves  in  nonlinear  WRF  simulations  ^  1),  we  assume 
that  nonlinear  dynamics  are  conhned  to  heights  z  <  Zi  nearer  the  mountain,  and  that  a 
nearly  linear  stationary  steady  wave  response  emerges  at  Zi  after  some  time  ti  (Lindeman 
et  ah  2008).  For  the  case  of  vertical  velocity,  WRF  output  initializes  the  FT  method  as 

w{k,  I,  Zi)  =  ww{k,  I,  Zi,  ti),  (8) 

where  ww{k,  I,  Zi,  ti)  is  the  horizontal  two-dimensional  Fourier  transform  of  the  WRF  vertical 
velocity  held  at  height  Zi  and  time  ti.  Final  values  of  Zi  and  ti  are  chosen  after  inspection 
of  the  WRF-simulated  dynamics  (see  section  2d):  typical  values  used  in  our  work  will  be 
Zi  =  6  km  and  ti  =  24  hours. 

The  gravity-wave  polarization  relations  allow  us  to  convert  w{k,  I,  z)  into  corresponding 
solutions  for  other  variables,  such  as  the  horizontal  velocity  components  u  and  v.  Here  we 
utilize  the  resulting  complex  spatial  solutions  u,  v,  and  w,  dehned  earlier,  whose  real  parts 
are  the  FT  solution  (4)  for  u,  v  and  w,  respectively,  and  whose  imaginary  parts  are  the 
Hilbert  transforms  of  those  FT  solutions,  thus  providing  local  amplitudes  and  phases  of  the 
velocity  solutions  (e.g.,  Mercier  et  ah  2008).  This  latter  property,  together  with  the  nvr  phase 
difference  between  horizontal  and  vertical  velocity  components  (n  G  Z) ,  allow  us  to  combine 
u,  V,  and  w  to  yield  local  solutions  of  phase-averaged  momentum  flux  per  unit  mass 


uw{x,  y,  z) 

1 

4 

{uw*  +  u*w) , 

(9) 

vw{x,  y,  z) 

1 

4 

{vw*  +  v*w) . 

(10) 

where  asterisks  denote  a  complex  conjugate.  The  relations  (9)  and  (10)  differ  from  other 
local  diagnostics  used  to  quantify  wave  momentum  fluxes  (e.g.,  Durran  1995)  by  using  a 


Fourier-ray  decomposition  of  the  WRF  output  to  synthesize  a  phase- averaged  result.  Their 
hdelity  depends  on  the  linear  stationary  wave  assumptions  underpinning  this  FT  method, 
and  so  are  assessed  when  applied  to  WRF  output  in  what  follows. 

c.  Drag 

1)  Surface  Pressure  Drag 

The  total  surface  pressure  drag  across  the  obstacle 

OO 

where  p'  is  the  surface  pressure  anomaly  (hereafter  we  denote  WRF  anomalies  with  lowercase 
primed  symbols  to  distinguish  them  from  spatial  realizations  of  the  Fourier-ray  solutions 
which  are  unprimed  lowercase). 

In  the  linear  {hm  -C  1)  hydrostatic  approximation,  (1)  and  (11)  yield  analytical  drag 
solutions  of  the  form 

|Dl|  =  Di(/?)  =  B{l5)apoNUhl,  (12) 

where  po  is  upstream  surface  density,  and  B{P)  is  the  scaling  form  factor  for  the  obstacle  (1) 
(Phillips  1984).  In  our  WRF  runs,  U  and  N  remain  constant  while  peak  mountain  heights 
hm  are  varied.  Thus,  the  linear  drag  (12)  scales  with  hm  =  Fr~^  in  our  WRF  runs  as 

IDJ  =  B(P)apo^hi,  (13) 

and  peaks  at  the  largest  (nonlinear)  hm  values. 

Consistent  with  previous  studies  (e.g.,  Miranda  and  James  1992;  Bauer  et  ah  2000),  it 
will  prove  most  useful  for  parameterization  applications  to  study  nondimensional  drag  values 
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normalized  by  the  relevant  linear  prediction:  e.g. 


Dp  =  —. 
^  Dr 


(14) 


In  compnting  Dp  =  |Dp|,  since  U  is  x-directed,  only  the  x-component  of  Dp  strictly  acts 
a  drag  on  the  flow  and,  due  to  the  symmetry  of  the  problem,  the  y  components  (sometimes 
referred  to  as  lift  forces)  should  integrate  to  zero.  Accordingly  many  studies  simply  compute 
the  x-component  of  (11)  only  (e.g.,  Epifanio  and  Durran  2001).  For  consistency,  here  we 
compute  both  components  and  evaluate  (11)  by  integrating  over  the  full  horizontal  WRF 
domain  of  — L/3  <  x  <  2L/3  and  — L/2  <  y  <  L/2. 


2)  Gravity- Wave  Contribution:  “Wave  Drag” 


The  component  of  the  drag  due  to  propagating  gravity  waves  at  a  given  height  z  is 


p-\-2L/3  p-\-L/2 

D^{z)  =  p{z)  /  /  {u'w' ,v'w')  dx  dy, 

J-L/3  J-L/2 


(15) 


where  u',  v'  and  w'  are  local  WRF  zonal,  meridional  and  vertical  velocity  anomalies,  respec¬ 


tively,  at  height  z.  In  the  linear  regime  one  should  recover  Dw{z)  =  Dp  k,  Dl  ai  all  heights, 
a  property  we  utilize  in  section  2d  to  tune  the  WRF  sponge  layer. 


The  FT  method  yields  corresponding  terms 


r-\-2L/3 

r+L/2 

Df,(s) 

=  p{z) 

l-L/3  . 

/-L/2 

{uw,  vw)  dx  dy, 

l*-\-2LfS 

r+L/2 

D,,(J) 

=  p(^) 

l-L/3  . 

' -L/2 

(uw,  vw)  dxdy. 

(16) 

(17) 


The  hrst  expression  is  the  FT  analogue  of  the  WRF  expression  (15)  using  individual  Fourier- 
ray  solutions  for  u,  v  and  w.  The  second  expression  integrates  the  phase-averaged  Fourier-ray 
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covariance  solutions  (9)  and  (10).  Since  (16)  and  (17)  give  identical  results,  for  simplicity 
here  we  only  use  Df^  and  refer  to  it  as  the  diagnostic  Fourier-ray  drag  Df.  We  also  compute 
a  linear  FT  drag  DfL  computed  from  (17)  using  the  linear  boundary  condition  (7). 

Since  these  expressions,  along  with  (11),  quantify  the  fraction  of  the  surface  pressure  drag 
Dp  that  appears  at  some  height  2;  as  gravity-wave  momentum  flux,  we  will  refer  to  them  as 
“wave  drag”  estimates.  This  is  not  to  be  confused  with  the  wave  drag  considered  in  OMD 
parameterizations,  which  is  the  decelerative  force  on  the  atmosphere  due  to  dissipation  of 
orographic  gravity-wave  momentum  flux  with  altitude,  which  we  do  not  study  in  this  paper. 


d.  Tuning  and  Validation 

A  series  of  linear  (h^  <C  1)  WRF  experiments  were  performed  with  varying  sponge  layer 
depths  and  damping  rates  until  computed  drag  reproduced  the  f  >  cx)  linear  prediction 
Dl  after  time  ti  without  subsidiary  structure  due  to  partial  reflections.  We  settled  upon  a 
vertical  damping  prohle  from  =15-30  km  based  on  eq.  (21)  of  Klemp  et  ah  (2008)  with  a 
damping  coefficient  at  z  =  30  km  of  0.2  s“^,  slightly  larger  than  the  value  chosen  by  Klemp 
et  ah  (2008)  due  to  our  deeper  sponge  layer  and  higher  horizontal  resolution. 

Figure  2a  plots  resulting  =  D^/ Dl  prohles  at  6  hourly  intervals  for  a  /3  =  1,  hm  =  0.1 
obstacle.  The  prohle  equilibrates  after  ~18  hours  to  values  close  to  the  hydrostatic  linear 
solution  Dl  =  1  with  little  evidence  of  spurious  structure  due  to  sponge-layer  rehections. 
Figure  2b  plots  the  corresponding  time  series  of  (normalized)  wave  drag  and  Df  at 
Zi  =  6  km,  along  with  surface  pressure  drag  Dp  and  the  linear  t  ^  00  predictions  D^  and 
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DfL.  While  surface  drag  Dp  — >  =  1  in  Figure  2b  after  only  a  few  hours,  wave  drag 

builds  up  over  time  to  a  limiting  value  of  ~  0.975  after  ~20  hours,  showing  that  ti  =  24  hours 
is  long  enough  to  generate  the  equilibrium  wave  response  at  z*  =  6  km. 

This  ~2%  reduction  in  from  in  Figure  2b  was  also  noted  by  Schar  and  Durran 
(1997)  in  their  hm  =0.05  simulation  for  an  axisymmetric  a  =  10  km  obstacle.  About  half  of 
this  difference  is  explained  by  the  ~1%  difference  between  the  linear  Fourier-ray  drag  Djl 
and  Dl  due  to  a  small  nonhydrostatic  reduction  in  drag  in  the  former  due  to  evanescent 
modes.  The  remaining  difference  is  due  to  lateral  propagation  of  wave  flux  through  the 
WRF  side  boundaries,  which  Djl  does  not  suffer  from  due  to  periodic  boundary  conditions 
and  wraparound.  This  boundary  absorption  increases  with  altitude  as  wave  helds  propagate 
upwards  and  disperse  more  horizontally,  as  evident  in  the  prohles  in  Figure  2a. 

Note  too  that  the  Df  and  curves  in  Figure  2b  lie  virtually  on  top  of  one  another, 
conhrming  the  hdelity  of  the  diagnostic  FT  relations  presented  in  section  2b  for  this  linear 
example.  Detailed  validation  of  these  FT  flux  diagnostics  for  other  linear  cases  simulated 
with  WRF  can  be  found  in  Lindeman  (2008). 


3.  Drag  Evolution  for  the  Axisymmetric  Obstacle 

Figure  3  plots  the  time  evolution  of  the  drag  due  to  the  axisymmetric  obstacle  (/3  =  1)  as 
hm  is  progressively  increased  from  0.1  to  1.5.  The  drag  values  are  normalized  by  the  linear 
hydrostatic  prediction  (12),  where  B{1)  =  n/A  (Phillips  1984). 

As  hm  — >■  1~,  the  equilibrium  surface  pressure  drag  Dp  becomes  progressively  larger  than 
unity  {Dl)-  In  the  hm  =0. 1-1.0  range  no  turbulent  kinetic  energy  (TKE)  is  generated  and 
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essentially  all  the  surface  drag  translates  into  propagating  gravity  wave  momentum  flux  at 
Zi  =  6  km  ^  Dp),  consistent  with  the  predictions  of  Figure  1.  However,  Figure  1  also 
shows  that  hm  values  slightly  larger  than  unity  enter  a  junction  point  where  all  the  dynamical 
regimes  meet  and  no  clear  dynamical  response  is  predicted. 

A  small  increase  in  hm  from  1.0  to  1.2  triggers  an  abrupt  change  in  simulated  drag 
behavior  in  Figure  3g.  Surface  drag  Dp  becomes  ~50%  larger  than  linear  predictions  and 
equilibrates  to  an  oscillating  state  with  a  period  of  close  to  4  hours.  Wave  drag  D^  and  Dj 
at  Zi  =  6  km  track  one  another  well  and  show  similar  4- hour  cycling  seen  in  Dp,  but  typically 
at  15-20%  lower  values  and  with  a  slight  phase  lag  due  to  vertical  propagation.  Thus  surface 
drag  no  longer  translates  fully  into  a  propagating  wave  response  at  6  km  altitude. 

The  periodic  nature  of  the  drag  is  due  to  cyclical  wave  buildup  and  breakdown,  as 
depicted  in  Figure  4  for  the  cycle  between  hours  20  and  24.  At  hour  20,  Dp  is  large  while 
Dy^  is  reduced  (Figure  3g)  as  forced  waves  propagate  upwards  and  build  in  amplitude  at 
Zi  =  6  km.  One  hour  later  wave  drag  Dy,  has  grown  to  the  point  where  isentropes  at  ~3  km 
in  Figure  4b  are  superadiabatic  and  generating  enhanced  TKE.  At  hour  22  wave  breaking 
peaks  with  a  large  burst  of  TKE  descending  the  obstacle  downslope.  This  TKE  has  migrated 
downstream  by  hour  23,  at  which  time  both  surface  drag  and  wave  drag  have  been  reduced 
by  accumulated  local  turbulent  dissipation  and/or  upstream  responses  to  wave  breaking  that 
increase  flow  splitting.  At  hour  24  the  isentropic  and  TKE  structures  in  Figure  4e  return  to 
a  morphology  almost  identical  to  that  at  hour  20  as  surface  drag  and  wave  activity  build  up 
again  and  the  cycle  repeats.  Note  that  TKE  (strong  nonlinearity)  is  conhned  below  Zi  =6  km, 
consistent  with  requirements  of  the  FT  method  and  with  the  agreement  between  Dy,  and 
Df  in  Figure  3e.  Despite  much  shorter  time  series,  axisymmetric  obstacle  simulations  of 
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Miranda  and  James  (1992)  and  Vosper  (2000)  show  evidence  of  similar  periodic  vacillations 
in  both  surface  and  wave  drag,  but  at  larger  hm  values. 

A  further  increase  of  to  1.5  yields  a  high-drag  state  with  wave  breaking  above  the 
lee  downslope.  In  contrast  to  our  =  1.2  result  and  the  =  1.5  results  of  Miranda  and 
James  (1992),  however,  drag  in  Figure  3h  does  not  vacillate  periodically.  Dp  instead  settles 
to  an  approximate  equilibrium  value  of  ~  1.3  while  Dy^  builds  up  gradually  to  an  equilibrium 
value  of  ~0.9.  The  former  compares  favorably  to  the  nonvacillating  equilibrium  Dp  ~1.2 
found  by  Schar  and  Durr  an  (1997). 

Figure  5  plots  local  gravity- wave  momentum  fluxes  at  Zi  =6  km  after  24  hours,  as  derived 
from  WRF  output  using  the  FT  relations  (9)  and  (10).  Figure  5a  reveals  the  expected  linear 
ship-wave  response  (Smith  1980)  at  hm  =  0.1,  with  net  westward  pmu{x,y,  z)  peaking  in 
magnitude  near  the  obstacle  peak,  and  smaller  pvw{x,  y,  z)  values  distributed  antisymmet¬ 
rically  about  y  =0  and  peaking  in  magnitude  in  the  divergent  wing  regions.  While  wave 
fluxes  at  h  =  1.0  in  Figure  5b  are  two  orders  of  magnitude  larger  [consistent  with  the  scaling 
(13)],  the  wave  pattern  is  very  similar,  except  the  region  of  maximum  p\uw\  is  now  displaced 
slightly  downstream.  At  hm  =  1-2  and  1.5,  wave  momentum  flux  patterns  are  noticeably 
modihed.  Peak  \pimv\  regions  are  narrower  and  located  above  the  obstacle  downslope.  For 
hm  =  1-5,  Lindeman  et  ah  (2008)  showed  that  this  results  from  leeside  collapse  and  ascent  of 
isentropes,  which  provide  a  sharp  deep  region  of  obstacle-like  forcing  to  the  flow  that  controls 
the  form  of  the  waves  at  6  km  more  than  the  actual  obstacle  shape  (see  also  Figure  4).  This 
narrower  forcing  generates  waves  with  shorter  horizontal  wavelengths  (Ding  et  ah  2003;  Lin¬ 
deman  et  ah  2008),  and  it’s  less  axisymmetric  nature  reduces  three-dimensional  propagation 
of  momentum  flux  into  diverging  wing  regions  in  Figures  5c  and  5d. 
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Figure  6  plots  the  drag  time  series  for  hm  =2-5.  Similar  responses  are  seen  in  all  cases. 
Dp  decreases  with  time  to  values  below  unity  and  does  not  equilibrate  after  24  hours.  By 
contrast,  steadily  increases  with  time  and  equilibrates  after  ~15  hours  at  values  signif¬ 
icantly  smaller  than  Dp,  consistent  with  greater  flow  splitting  and  lee-vortex  generation  at 
these  larger  hm  values  (Schar  and  Durran  1997). 

A  cross  section  of  isentropes  and  TKE  at  j/  =  0  after  24  hours  is  shown  in  Figure  7 
for  the  hm  =  3  simulation.  A  sharp  zone  of  isentropic  collapse  and  ascent  is  evident  above 
the  top  third  of  the  lee  downslope,  with  TKE  peaked  sharply  within  this  leeside  trough 
at  unstable  temperature  gradients.  As  for  the  hm  =1.5  obstacle  (Lindeman  et  al.  2008), 
Figure  7  shows  that  this  leeside  isentropic  structure  provides  the  major  obstacle-like  forcing 
for  the  waves  that  emerge  a.t  Zi  =  6  km,  which  are  concentrated  downstream  of  the  obstacle 
peak.  Similar  leeside  offsets  in  wave  activity  are  evident  in  the  hm  =  3.0  results  of  Schar  and 
Durran  (1997).  Since  the  leeside  isentropic  structure  in  Figure  7  is  quite  steady  with  time 
(not  shown),  the  wave  fields  it  forces  are  approximately  stationary  and  steady,  and  thus  their 
integrated  momentum  fluxes  D^  are  well  described  by  the  FT  diagnostic  Dj  in  Figure  6b. 

Figure  8  shows  the  corresponding  local  wave  momentum  fluxes  puw  and  pvw  at  6  km. 
Zonal  fluxes  peak  in  magnitude  above  the  obstacle  downslope  due  to  stationary  obstacle-like 
forcing  by  lee-side  isentropic  structures  as  in  Figure  7.  These  peak  fluxes  are  concentrated 
within  regions  of  progressively  smaller  horizontal  extent  as  hm  increases,  consistent  with 
a  progressively  higher  dividing  streamline  and  flow  above  passing  over  a  smaller  narrower 
uppermost  obstacle  peak,  forcing  waves  of  correspondingly  smaller  horizontal  wavelength  (see 
Figure  7).  At  the  largest  (most  nonlinear)  hm  =  ^  and  5,  subsidiary  \pWw\  maxima  appear 
in  Figure  8c  and  8d  at  the  lateral  edges  of  the  obstacle,  due  to  flow  splitting  and  diversion 
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around  the  obstacle  peak,  resulting  in  Bernoulli-type  accelerations  that  force  waves  at  the 
cross-stream  extremities  of  the  obstacle.  A  similar  high-hm  transition  in  wave  forcing  to  the 
onter  edges  of  the  obstacle  was  seen  in  the  (3  =  5  simnlations  of  Olafsson  and  Bongeanlt 
(1996).  The  fine-scale  wave  strnctnres  at  the  largest  hm  may  be  compromised  somewhat 
by  onr  1.5  km  gridpoint  resolntion  as  the  clipped  obstacle  width  for  wave  forcing  becomes 
progressively  narrower  as  the  dividing  streamline  gets  progressively  higher  with  increasing 
hm- 

Since  some  drag  responses  do  not  obvionsly  eqnilibrate  after  24  hours,  the  hm  =  1-2  and 
3.0  WRF  simnlations  were  integrated  ont  to  48  honrs.  Figure  9  shows  the  resulting  drag 
prohles.  Oscillatory  pressnre  and  wave  drag  at  hm  =  1.2  is  clearly  a  robnst  eqnilibrinm 
resnlt.  The  secnlar  decrease  in  Dp  ont  to  24  honrs  seen  for  hm  =  3.0  abates  in  Fignre  9b 
after  ~28  honrs,  wherepon  Dp  increases  slightly  then  eqnilibrates  at  ~0.5  at  ~36  honrs. 
Wave  flnxes  remain  largely  nnchanged  dnring  this  period.  These  integrations  snggest  that 
drag  in  the  20-24  honr  range  is  close  to  the  long-time-limit  eqnilibrinm  response. 

Thns  Fignre  10  plots  mean  drag  as  a  fnnction  of  hm  over  the  20-24  honr  eqnilibrinm 
range,  both  in  normalized  and  absolnte  form.  At  0  <  hm  <  1,  Fignre  10a  shows  a  gradnal 
linear  increase  in  Dp  above  nnity.  A  similar  small  linear  increase  over  this  hm  range  was 
modeled  by  Epifanio  and  Dnrran  (2001)  using  weakly-nonlinear  semianalytical  theory.  This 
linear  trend  breaks  at  hm  ^  1  with  the  onset  of  a  peak  oscillating  high-drag  wavebreaking 
state  of  Dp  =1.4-1.55  at  hm  =  1.2,  with  Dp  decreasing  monotonically  with  hm  thereafter. 
This  mean  Dp  variation  with  hm  agrees  with  that  fonnd  by  Baner  et  ah  (2000)  (who  nsed 
a  slightly  different  obstacle  shape)  bnt  differs  somewhat  from  the  resnlts  of  Miranda  and 
James  (1992),  who  fonnd  Dp  peaking  at  ~1.8  at  hm  ~  1.5.  Onr  Dp  valnes  at  hm  =1.5 
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and  3  are  similar  to  those  of  Schar  and  Dnrran  (1997)  bnt  are  again  smaller  than  those  of 
Miranda  and  James  (1992).  Onr  resnlts  snpport  the  contention  of  Baner  et  ah  (2000)  that 
the  Miranda  and  James  (1992)  drag  valnes  are  systematically  too  large. 

Next  we  stndy  mean  valnes  and  how  they  relate  to  snrface  pressnre  drag  in  the 
nonlinear  hm  range  where  D^/Dp  <  1.  A  common  parameterization  approach  is  to  estimate 
the  height  Zds  of  the  dividing  streamline,  separating  flow  over  the  obstacle  peak  {z  >  Zds) 
from  flow  blocked  or  diverted  aronnd  the  obstacle  (z  <  Zds),  as 

Zds  =  (1  -  aFr)  hm  =  (^hm  -  (18) 

(Hnnt  and  Snyder  1980),  where 

a  =  hcs  =  Fr~^,  (19) 

the  critical  inverse  Fronde  nnmber  for  flow  splitting.  The  “clipped”  obstacle  for  wave  forcing 
is  then  hc{x,y)  =  ma.x[h{x ,  y) ,  Zds]  (Lindeman  et  ah  2008)  which  has  a  peak  height 

^ds 


whose  normalized  height,  from  (2),  is  eqnal  to  a.  Fignre  1  predicts  that  a  varies  with  /3, 
and  it  can  also  vary  with  different  dehnitions  of  Zds  (see  section  6.6  of  Baines  1995).  Here, 
for  simplicity,  we  will  adopt  the  Sheppard  criterion  of  a  =  1  which,  despite  its  theoretical 
shortcomings  (Smith  1989b),  reprodnces  laboratory  measnrements  of  Zds  (Snyder  et  ah  1985; 
Vosper  et  al.  1999;  Ding  et  ah  2003). 

Large  triangles  in  Fignre  10a  show  the  resnlting  wave  drag  D^.  calcnlated  by  the  FT 
method  nsing  hc{x,y)  as  the  linear  bonndary  condition  (7)  at  Zi  =  Zds-  An  even  simpler 
analytical  approximation  assnmes  the  linear  solntion  (12)  for  the  clipped  obstacle  height 
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{hc)m,  which  yields  a  constant  “saturated”  clipped  wave  drag  (see  Appendix  A)  of 


Del  =  B{P)apoa‘ 


N' 


(21) 


plotted  with  small  triangles  in  Figure  10b. 

Both  Dc  and  Dd  systematically  underestimate  the  actual  wave  drag  at  hm  >  1  in  Fig¬ 
ure  10.  Nonetheless,  Dd  reproduces  the  flat  trend  in  D^  at  hm  >  2.5  in  Figure  10b  but  with 
Dyj/ Dd  ~  1.5:  an  a  of  1.2  instead  of  1  would  make  this  ratio  near  unity,  a  choice  that  can 
be  justihed  from  Figures  1,  3  and  10a.  D^/Dd  is  larger  in  the  range  hm  =1.2-2  with  a  peak 
near  2,  demarking  this  region  as  a  “high  wave-drag  state.” 

It  is  also  useful  to  study  how  D^  might  scale  with  Dp  given  that  many  OMD  param- 
eterizations  quantify  gravity-wave  fluxes  by  scaling  surface  drag.  A  reasonable  Up-scaled 
empirical  fit  to  D^  at  all  hm,  shown  as  the  solid  gray  curve  in  Figure  10,  is: 


D^  ~  D fit 


Dp 

0  <  hm  <  1, 

Dphm^-^ 

1  <  hm  <  2.5, 

hm  >  2.5. 

(22) 


4.  Drag  Evolution  for  the  Elliptical  Obstacle 

Figure  11  plots  time  evolution  of  the  elliptical  {(5  =  3)  obstacle  drag  as  hm  is  increased 
from  0.1  to  1.5.  Again,  drag  is  normalized  by  the  corresponding  hydrostatic  linear  prediction 
Dl  in  (12),  where  5(3)  ~  0.9369/9  ~  2.82  (Phillips  1984).  Thus  the  linear  drag  in  this  case 
is  enhanced  over  that  of  the  axisymmetric  obstacle  by  Di[?))/Dl{1)  =  5(3)/5(l)  3.6. 

At  /3  =  3,  Figure  1  predicts  nonbreaking  waves  at  0  <  hm  <  heb  ~  0.73,  followed  by  an 
extended  wave  breaking  regime  at  hm  >  heb-  Drag  at  hm  =0. 1-0.5  in  Figure  11  reveals  pure 
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wave  responses,  for  which  Dp  Ri  D^j  ^  Dj  Di,  that  generate  no  TKE.  The  =  0.7  drag 
in  Fignre  lid  oscillates  dne  to  small  vacillating  levels  of  TKE  (peaking  at  ~0. 01-0. 05  s“^) 

that  form  at  nnstable  wave  phases,  consistent  with  the  predicted  wave  breaking  threshold 
of  hch  ~  0.73  in  Fignre  1.  At  hm  =0.85-1.5,  strong  wave  breaking  dynamics  are  simnlated, 
again  consistent  with  Fignre  1.  The  snrface  drag  Dp  at  hm  =0.85-1.5  nndergoes  large  cyclical 
excnrsions  with  a  period  that  decreases  from  ~9  honrs  at  hm  =  0.85  (Fignre  lie)  to  ~6  honrs 
at  hm  =  1-5  (Fignre  llh). 

Unlike  the  corresponding  (3=1  resnlts  in  Fignre  3,  the  diagnostic  FT  wave  drag  Df 
does  not  reprodnce  Dyj  well  at  hm  ^  0.7  in  Fignre  11.  Reasons  why  are  highlighted  in 
Fignre  12,  which  plots  isentropes  and  TKE  cross  sections  at  2  honrly  intervals  along  the 
y  =  0  axis  for  the  hm  =  0.85  simnlation.  At  honr  3  when  Df  and  agree  in  Fignre  lie, 
nonlinearity  in  the  form  of  overtnrning  isentropes  and  large  levels  of  TKE  is  conhned  below 
Zi  =  6  km  (where  Df  and  D^  are  evalnated),  consistent  with  the  assnmptions  of  the  FT 
method  (see  section  2b).  However,  nnlike  the  (3  =  1  cases  (see  Fignre  4),  nonlinearity  is  not 
conhned  below  6  km  at  later  times.  For  example,  at  honr  7  there  is  strong  overtnrning  and 
tnrbnlent  mixing  at  ~10  km  altitnde  (Fignre  12c).  At  later  times  high-altitnde  wave  activity 
is  evident  well  downstream  of  the  obstacle  in  Fignres  12d  and  12e.  The  sonrces  of  these  waves 
are  near-snrface  wake  dynamics  (e.g.,  Dnpont  et  al.  2001)  and  wavebreaking  regions  at  higher 
altitndes  (e.g.,  Bacmeister  and  Schoeberl  1989;  Satomnra  and  Sato  1999;  Tan  and  Eckermann 
2000).  The  resnlting  waves  are  generally  nonstationary  and  some  may  propagate  downwards 
(Bacmeister  and  Schoeberl  1989),  and  thns  are  not  properly  characterized  by  the  stationary 
npgoing  gravity-wave  dispersion  relation  (6)  that  nnderpins  the  FT  diagnostics. 

Drag  at  the  larger  hm  valnes  of  2-5  is  plotted  in  Fignre  13.  The  vacillating  wave  breaking 
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regime  extends  to  at  least  hm  =  3,  where  the  drag  vacillation  is  now  weak  and  its  period 
has  rednced  to  ~4  honrs  (Fignre  13b).  At  still  larger  hm  the  pressnre  drag  decreases  more 
nniformly  with  time  and  wave  drag  settles  to  eqnilibrinm  values,  similar  to  the  behavior 
found  for  the  /3  =  1  obstacle  at  similarly  large  hm  values  in  Figure  6. 

Given  the  large  drag  vacillations  for  the  (3  =  3  obstacle,  mean  drag  is  computed  over  the 
last  12  hours  of  each  simulation,  with  resulting  means  and  standard  deviations  plotted  in 
Figure  14.  Mean  drag  in  the  range  0  <  hm  ^  0.5  does  not  exhibit  a  noticeable  linear  increase 
with  hm  as  was  found  for  the  axisymmetric  obstacle,  but  instead  tracks  the  hydrostatic 
analytical  value  of  Dl  =  1  closely.  This  persists  to  hm  =  0.7  despite  the  vacillation  in 
Figure  lid  and  associated  standard  deviation  in  Figure  14.  At  hm  =  0.85,  a  large  abrupt 
increase  occurs  in  both  Dp  and  its  standard  deviation,  due  to  sudden  onset  of  a  vacillating 
high-drag  wavebreaking  state.  The  high-drag  Dp  enhancements  peak  at  hm  ~  1.  At  larger 
hm,  mean  Dp  decreases  monotonically  as  do  the  associated  standard  deviations.  Wells  et  ah 
(2008)  report  a  mean  Dp  ~  1.1  at  hm  =  2  in  their  Boussinesq  simulations  without  Coriolis 
terms,  somewhat  larger  than  our  mean  value  here  of  0.83  ±  0.13. 

Wave  drag  .D^,  which  tracks  Dp  and  Dl  at  small  hm,  drops  abruptly  below  Dp  at  hm  ^ 
0.85  as  wave  breaking  and  (at  larger  hm)  how  splitting  cause  progressively  greater  fractions 
of  the  surface  pressure  drag  to  manifest  as  near-surface  dynamics  rather  than  gravity  waves. 
Despite  the  strongly  nonlinear  dynamics  that  yield  a  generally  nonstationary  nonlinear  wave 
held  at  Zi  =  6  km  (Figure  12),  the  unnormalized  areally-averaged  momentum  hux  is 
htted  fairly  well  throughout  the  extended  nonlinear  range  hm  ~0.85-5  in  Figure  14b  by  the 
simple  clipped  obstacle  prediction  D^i  in  (21).  For  example,  mean  D^  values  at  hm  =3-5  are 
only  ~10-20%  larger  than  D^.  At  hm  =0.85-3,  D^  exceeds  D^i  by  more  than  this  factor. 
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consistent  with  this  range  of  hm  being  a  high-drag  wave-breaking  state. 


5.  Discussion 

Figure  15  compares  the  mean  normalized  surface-pressure  and  far-held  wave  drags  for 
the  /9  =  1  and  (3  =  3  obstacles  as  a  function  of  hm-  While  =  Dp  for  both  obstacles  in  the 
mutual  “linear”  (non-TKE)  range  hm,  ~0.1-0.7,  the  Dy^  curve  for  the  axisymmetric  obstacle 
increases  linearly  with  hm  in  this  range  to  values  progressively  larger  than  unity,  whereas 
the  Dw  curve  for  the  elliptical  (3  =  3  obstacle  stays  close  to  the  hydrostatic  linear  prediction 
of  Dl  =  1. 

In  the  high-drag  region  of  hm  ~l-2,  mean  Dp  values  for  the  (3  =  3  obstacle  peak  at  a 
smaller  hm  and  attain  larger  values  than  for  the  axisymmetric  obstacle,  and  individual  values 
exhibit  a  much  larger  variation  with  time  (standard  deviation).  This  is  broadly  consistent 
both  with  Figure  1  and  the  more  intense  wave  breaking  (relative  to  axisymmetric  obstacles) 
expected  for  elliptical  obstacles  (1)  whose  long  axis  b  is  aligned  perpendicular  to  the  flow, 
due  to  less  dispersion  of  wave  activity  that  leads  to  more  concentrated  and  unstable  wave 
helds.  The  smaller  D^  values  at  =  6  km  for  the  (3  =  3  compared  io  (3  =  1  obstacle 
are  consistent  with  greater  near-surface  dissipation  of  waves  by  overturning  and  TKE  in  the 
(3  =  3  simulations  (c.f.  Figures  4  and  12). 

In  the  nominal  flow-splitting  regime  of  hm  ^  2,  Dp  values  in  Figure  15  are  almost 
identical  for  both  obstacles.  Using  a  different  obstacle  shape,  Bauer  et  ah  (2000)  also  found 
very  similar  Dp  values  at  large  hm  for  (3  =1-4.  Some  analytical  models  of  flow-blocking  drag 
predict  Dp  oc  hm  at  large  hm  (Lott  and  Miller  1997;  Scinocca  and  McFarlane  2000),  while 
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others  suggest  Dp  Di  ~  1  due  to  surface  drag  enhancements  from  rotation  and  friction 
(Webster  et  ah  2003).  Our  common  Dp  curve  is  well  £t  as  at  hm  =2-5.  Stein  (1992) 
reported  a  Dp  £t  oc  in  this  range  based  on  two-dimensional  simulations. 

Absolute  wave  drag  D^  at  hm  ~2.5-5  was  almost  constant  for  both  the  (3  =  1  and  (3  =  3 
obstacles,  and  predicted  reasonably  accurately  by  scaling  the  “saturated”  linear  wave  drag 
Del  given  by  (21)  (see  Figures  10b  and  14b).  It  is  remarkable  that  Dy^  values  are  all  so 
similar  in  this  range,  given  the  large  variation  in  local  wave  momentum-flux  morphology  as 
a  function  of  hm  in  Figure  8,  and  nonstationary  transient  waves  from  surface  wake  and  TKE 
sources  for  the  (3  =  3  obstacle  in  Figure  12.  While  this  ~  D^  hnding  may  prove  to  be 
very  convenient  for  parameterizing  orographic  gravity  wave  momentum  fluxes  in  this  range, 
it  must  be  recognized  that  the  simple  theoretical  arguments  that  yield  Dd  (see  Appendix 
A)  are  inconsistent  with  the  WRF-simulated  wave  helds.  As  explicitly  demonstrated  by 
Lindeman  et  al.  (2008),  the  effective  obstacle  forcing  for  these  waves  does  not  resemble  the 
clipped  obstacle  forcing  underpinning  Dd-  It  is  instead  dominated  by  obstacle-like  forcing 
produced  by  the  leeside  fluid  dynamical  response,  specihcally  the  sharp  descent  and  ascent  of 
isentropes  along  the  obstacle  downslope  (see  Figure  7),  which  yield  highly  concentrated  wave 
momentum  flux  in  Figure  8  above  the  obstacle  downslope,  rather  than  near  the  peak.  In 
short,  (21)  gives  an  approximately  correct  answer  but  for  demonstrably  wrong  reasons.  Why 
Dd  is  a  reasonable  £t  to  modeled  wave  fluxes  in  this  nonlinear  range  is  not  yet  clear,  but  this 
need  not  be  a  barrier  to  using  it  in  parameterizations.  For  example,  parameterizations  based 
the  saturated  gravity-wave  spectra  observed  in  the  middle  atmosphere  are  now  common  in 
models,  despite  the  lack  of  an  accepted  theory  for  how  these  spectra  are  formed  (McLandress 
and  Scinocca  2005). 
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Indeed,  several  OMD  parameterizations  already  use  relations  akin  to  Dd  to  assign  wave 
momentum  fluxes  in  this  nonlinear  range  (Scinocca  and  McFarlane  2000;  Webster  et  al. 
2003).  The  value  of  a,  defined  equivalently  in  these  OMD  schemes  via  (19)  as  the  critical  in¬ 
verse  Fronde  number  Fr~^  (since  Fr~^  ^  is  also  assumed),  must  be  chosen.  Originally 
these  schemes  used  Fr~^  ~0.5-l  (Scinocca  and  McFarlane  2000;  Webster  et  al.  2003)  but, 
given  a  lack  of  theoretical  closure  on  D^/Dp  values  at  hm  >  Fr~^,  Fr~^  has  been  used  as  a 
“free”  tuning  parameter  in  models  and  subsequently  revised  downwards  to  ~0. 25-0. 375  to 
improve  skill  (Wells  et  al.  2005;  Scinocca  et  al.  2008).  Given  our  finding  of  Dw/Dd  ~1. 1-2.0 
for  a  Dd  calculation  using  =  a  =  1  and  the  Fr~^  dependence  of  Dd  in  (21),  these  re¬ 
sults  imply  order  of  magnitude  underestimates  of  orographic  gravity  wave  momentum  fluxes 
at  hm  ^  1  for  Ftc  ~0. 25-0. 375. 

This  is  illustrated  further  in  Figure  16,  which  plots  D^j/Dp  as  a  function  of  hm  for  both 
the  P  =  1  and  P  =  3  obstacles.  Some  previous  simulations  have  suggested  that  D^/Dp 
rapidly  decreases  at  hm  ^  Fr~^  and  essentially  vanishes  at  large  hm  (Miranda  and  James 
1992;  Wells  et  al.  2005).  While  Dyj/Dp  certainly  drops  suddenly  below  unity  at  hm  ~  Fr~^ 
in  Figure  16,  especially  for  the  P  =  3  obstacle,  our  results  also  show  that,  even  at  large 
hm  ~2-5,  propagating  “far  held”  gravity  waves  remain  a  significant  (>20%)  fraction  of  the 
total  pressure  drag.  Our  hm  =  2,  /3  =  3  value  of  0.43  ±  0.10  at  Zj  =  6  km  in  Figure  16b 
compares  reasonably  to  a  corresponding  value  of  ~0.52  evaluated  by  Wells  et  al.  (2008)  at 
~100  m  above  their  obstacle  peak.  The  slightly  smaller  mean  value  here  is  consistent  with 
some  additional  near-held  turbulent  dissipation  of  waves  in  propagating  to  =  6  km. 

Between  these  asymptotic  wave-flux  limits  of  Dyj  ^  at  hm  <C  1  and  ~  Dd 
at  hm  3>  1,  Figures  10b  and  14b  also  show  an  intermediate  zone  at  hm  ~l-2.5  where 
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attains  maximum  mean  values  and  can  vacillate.  By  analogy  to  surface  drag,  we  refer  to  this 
region  as  a  “high- wave- drag”  state.  However,  unlike  surface  drag,  does  not  exceed  in 
this  high-drag  range,  but  instead  exceeds  the  saturated  wave-drag  Dd,  or  more  accurately, 
exceeds  the  limiting  at  3>  1.  These  wave-drag  enhancements  are  fairly  small  for  the 
P  =  3  obstacle  (presumably  due  to  greater  low-level  dissipation  due  to  wave  breaking)  but 
are  on  the  order  of  50%  for  the  axisymmetric  obstacle  in  Figure  10b. 

For  OMD  parameterization,  our  isolation  of  the  far-held  gravity-wave  component  of 
the  total  surface  drag  Dp  in  turn  dehnes  the  component  Db  =  Dp  —  D^  that  is  conhned  and 
deposited  below  Zi  =  6  km.  Figure  17  plots  its  variation  with  hm  in  both  normalized  and 
absolute  form  for  both  the  /3  =  1  and  (3  =  3  obstacles.  While  Db  is  usually  referred  to  as 
how-blocking  drag  in  OMD  parameterizations,  here  it  likely  includes  other  ehects,  such  as 
hux  deposition  due  to  low-level  wave  breaking. 

Since  this  study  is  a  hrst  methodical  attempt  to  characterize  far-held  gravity-wave  com¬ 
ponents  of  the  drag  from  three-dimensional  obstacles  at  a  range  of  nonlinear  hm  values,  we 
employed  a  number  of  simplihcations  to  make  wave  characterization  easier.  One  potentially 
important  omission  was  the  Coriolis  force  (/  =  0).  Olafsson  and  Bougeault  (1997)  compared 
surface  drag  in  /  =  0  and  /  %  0  simulations  and  found  larger  Dp  values  at  large  hm  for 
/  %  0.  They  used  a  =  40  km  rather  than  the  10  km  value  used  here.  Wells  et  ah  (2008) 
argued  that  the  primary  ehect  of  nonzero  /  at  obstacle  widths  0(10  km)  was  to  introduce 
asymmetries  to  the  how,  which  trigger  vortex  shedding  that  yields  larger  Dp  values  (e.g. 
Vosper  2000).  For  example.  Wells  et  ah  (2008)  found  similar  Dp  between  /  %  0  simulations 
and  /  =  0  simulations  in  which  asymmetry  was  introduced  artihcally  with  a  potential  tem¬ 
perature  perturbation  (Schar  and  Durran  1997).  We  included  neither  rotation  nor  potential 
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temperature  perturbations  in  these  runs  to  keep  our  study  of  wave  momentum  fluxes  as  sim¬ 
ple  as  possible.  As  a  result,  our  drag  values  may  be  underestimated  at  large  hm-  That  having 
been  said,  our  nonlinear  simulations  yield  asymmetric  lee  vortex  evolution  after  ~18  hours. 

For  the  /3  =  3  obstacle,  an  additional  question  for  OMD  parameterization  is  how  wave 
momentum  fluxes  at  hm  ^  Fr~^  vary  as  the  angle  (p  between  the  incident  flow  vector  and 
the  short  axis  of  the  obstacle  changes  from  our  value  of  =  0°.  While  unable  to  predict 
mean  magnitudes.  Wells  et  ah  (2008)  found  that  the  linear  drag  relations  of  Phillips  (1984) 
provided  a  good  £t  to  how  surface  pressure  drag  varied  as  a  function  of  (p  for  a  /3  =  3  obstacle 
at  hm  =  2.  They  also  found  a  strong  D^/ Dp  variation  with  (p,  from  ~0.5  at  (p  =  0°  to  values 
~0.9  at  if  near  ±90°.  Since  the  latter  p  values  are  equivalent  to  changing  from  a  /3  =  3 
to  a  =  1/3  obstacle.  Figure  1  predicts  for  hm  =  2  and  (3  =  1/3  that  the  flow  is  near  the 
threshold  between  linear  wave  and  flow-splitting  dynamics,  for  which  D^/Dp  ~  1  would  be 
expected  based  on  the  current  results. 

Other  omitted  processes  may  also  be  important  in  setting  gravity  wave  momentum  fluxes. 
For  example,  using  a  no-slip  lower  boundary  condition,  Olafsson  and  Bougeault  (1997)  found 
surface  drag  reductions  at  hm  ~l-3  and  enhancements  at  larger  hmi  relative  to  correspond¬ 
ing  free-slip  results,  although  subsequent  simulations  suggest  smaller  frictional  modihcations 
at  large  hm  (Wells  et  ah  2005).  From  our  OMD  parameterization  perspective,  we  can  justify 
the  omission  of  surface  friction  here  by  noting  that  it  is  parameterized  separately  in  global 
models  with  orographic  roughness,  land-surface  and  planetary  boundary  layer  parameteri- 
zations.  Thus,  frictional  modihcations  to  parameterized  OMD  would  ideally  involve  inputs 
from  those  schemes  that  modify  a  frictionless  baseline  OMD  parameterization.  In  practice, 
however,  current  OMD  parameterizations  in  models  are  not  modihed  by  frictional  inputs. 
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As  a  result,  the  model-tuned  values  of  Fr~^  ~0. 25-0. 375  used  in  some  OMD  parameteriza- 
tions,  which  signihcantly  underestimate  the  gravity  wave  momentum  fluxes  reported  here, 
may  parameterize  the  mean  overall  effect  of  friction  and  other  omitted  processes  in  reducing 
far-held  wave  huxes.  Thus,  future  investigations  of  the  inhuence  of  friction,  (p,  and  other 
omitted  processes  on  wave  momentum  huxes  at  nonlinear  Froude  numbers  would  clearly  be 
worthwhile. 


6.  Summary 

Using  WRF  model  simulations,  we  studied  mesoscale  dynamics  resulting  from  a  steady 
unsheared  how  impinging  upon  an  idealized  obstacle  (1)  of  along-stream  width  a  =  10  km 
and  cross-stream  width  b  =  f3a.  We  focused  on  an  axisymmetric  {f3  =  1)  and  elliptical 
(/3  =  3)  obstacle,  and  for  each  performed  experiments  in  which  the  normalized  obstacle 
height  hm  =  Fr~^  was  progressively  increased  from  the  linear  regime  {hm  -C  1)  through 
to  the  highly  nonlinear  regime  {hm  3>  1),  as  shown  in  Figure  1.  The  primary  foci  of  our 
analysis  were  the  local  and  areally-integrated  momentum  huxes  of  the  far-held  equilibrium 
gravity  wave  helds  that  propagated  to  a  height  of  6  km.  Those  huxes  were  quantihed  both 
directly  from  the  WRF  velocity  anomalies  and  indirectly  using  diagnostic  Fourier-transform 
methods.  We  sought  mean  results  that  could  be  used  to  specify  gravity-wave  momentum 
huxes  in  subgrid-scale  orographic  mesoscale  drag  (OMD)  parameterizations  in  global  weather 
and  climate  models.  Our  major  hndings  are  summarized  schematically  in  Figure  18. 

Since  many  OMD  parameterizations  assign  wave  momentum  huxes  by  scaling  surface 
pressure  drag.  Figure  18a  hrst  summarizes  the  bulk  properties  of  (normalized)  pressure  drag 
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as  a  function  of  h^.  At  /im  <  Fr^  ^  pressure  drag  is  close  to  linear  hydrostatic  predictions 
(Dr  ~  1)  and  translates  fnlly  into  propagating  gravity  wave  momentnm  flnx  at  6  km.  At 
Fr~^  <  hm  ^  2  high-drag  wave-breaking  states  {Dp  >  1)  occnr,  with  the  Fr~^  threshold 
predicted  to  within  ~10%  for  both  the  /3  =  1  and  /3  =  3  obstacle  nsing  the  Smith  (1989b) 
cnrves  in  Fignre  1.  Snrface  drag  for  the  (3  =  3  obstacle  vacillates  in  the  hm  ~0.7-3  range 
dne  to  cyclical  bnildnp  and  breakdown  of  wave  activity,  with  the  period  and  amplitnde  of 
these  vacillations  decreasing  with  increasing  hm-  Periodic  drag  vacillation  occnrred  for  the 
axisymmetric  obstacle  only  at  hm  =  1-2.  At  hm  ^2-3  pressnre  drag  decreased  smoothly 
with  time  to  normalized  eqnilibrinm  valnes  that  were  almost  identical  for  both  obstacles, 
and  were  well  htted  as 

Mean  featnres  in  the  areally-integrated  wave  momentnm  flnxes,  or  “wave  drag,”  at  Zi  = 
6  km  are  snmmarized  in  absolnte  (nnnormalized)  form  in  Fignre  18b.  In  the  linear  hm  <  Fr~^ 
range,  wave  drag  matches  both  Dp  and  Dl  and  thus  varies  according  to  the  latter’s 
analytical  dependence  in  (13).  While  at  hm  >  Fr~^  suddenly  becomes  a  signihcantly 
smaller  fraction  of  the  total  pressure  drag  (see  Figure  16),  absolute  values  peak  in  the  high- 
drag  range  Fr~^  <  hm  ^2.  As  with  pressure  drag,  wave  drag  vacillates  in  this  range  due  to 
periodic  buildup  and  breakdown  of  wave  amplitudes  just  above  the  obstacle  (see  Figure  4). 

At  hm  ^2-3  wave  drag  levels  off  to  a  constant  value,  which  is  quantihed  to  within  10- 
40%  by  a  “saturated”  clipped  wave  drag  Dd  given  by  the  analytical  relation  (21)  using 
a  =  Fr~^  K.  1.  Stable  invariant  wave  drag  in  this  high-hm  range  is  a  surprising  hnding, 
given  that  the  three-dimensional  distribution  of  wave  momentum  flux  varies  appreciably  with 
hm  over  this  range  (see  Figure  8)  and  can  be  additionally  impacted  by  nonstationary  wave 
responses  in  the  elliptical  obstacle  simulations  (see  Figure  12).  We  show  that  the  analytical 
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relation  Dd  arises  here  from  a  demonstrably  incorrect  model  of  the  h^,  ^  Fr~^  wave  held. 
Nonetheless,  based  on  the  agreement  in  Fignre  18b,  it  offers  a  simple  and  convenient  way 
of  accurately  quantifying  wave  hux  in  this  nonlinear  hm  range.  Mean  wave-drag  curves  like 
those  in  Figure  18b  can  easily  be  incorporated  into  existing  OMD  parameterizations,  given 
their  close  current  associations  with  the  mean  drag  properties  of  idealized  three-dimensional 
obstacles  (see,  e.g.,  Scinocca  and  McFarlane  2000). 

Indeed,  some  current  OMD  parameterizations  already  use  analogues  of  ~  Dd  to 
assign  wave  momentum  huxes  at  nonlinear  Fronde  numbers.  While  we  found  best  agreement 
using  Fr~g^  =1-1.2,  the  values  of  Fr~^  Fr~J^)  used  in  OMD  parameterizations  that 
have  been  “tuned”  to  optimize  model  skill  are  currently  ~0. 25-0. 375,  which,  based  on  our 
results,  imply  order-of-magnitude  underestimates  of  parameterized  orographic  gravity  wave 
momentum  huxes  at  hm  >  Fr~^  in  such  models.  However,  our  model  simulations  do  not 
include  additional  ehects  like  surface  friction  and  rotation  which  may  reduce  wave  huxes. 
Our  results  also  show  large  vacillations  in  both  surface  and  wave  drag,  particularly  in  the 
range  Fr~^  <  hm  ^2-3,  which  existing  schemes  do  not  include,  but  could,  for  example,  via 
an  explicitly  stochastic  parameterization  of  OMD  (see,  e.g..  Palmer  2001). 
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APPENDIX  A 


Analytical  Wave  Drag  Relations  Based  on  the  Dividing 

Streamline  Concept 

We  consider  waves  forced  by  flow  over  a  “clipped”  obstacle  hc{x,y)  located  between  the 
dividing  streamline  Zds  and  the  obstacle  peak  hm  (Lindeman  et  ah  2008).  The  simplest 
derivation  assnmes  the  same  shape  as  the  nnclipped  obstacle,  same  widths  a  and  b,  bnt 
reduced  height  {hc)m  =  hm  —  Zds  given  by  (20).  Assnming  the  linear  hydrostatic  solntion  Dl 
holds  at  0  <  hm  <  o,  then  aX.  hm>  a  the  wave  drag  is 


(Al) 


In  reality  the  clipped  monntain  peak  is  narrower  and  generates  mnch  shorter-scale  waves 
(Ding  et  al.  2003;  Lindeman  et  al.  2008).  A  better  approximation  then  might  be  to  assnme 
that  the  effective  half  widths  Oc  and  be  of  the  clipped  obstacle  hc{x,y)  rednee  linearly  with 
the  clipping  height  so  that  the  obstacle’s  height-to- width  ratio  remains  constant.  Ignoring 
possible  nonhydrostatic  effects  as  Oc  and  be  contract,  this  yields  a  clipped  wave-drag 


Do2  =  Di,a:'hJ. 


(A2) 


In  fact,  neither  approximation  hts  the  actnal  solntion  De  from  the  Fonrier-ray  calcnlation 


(large  triangles  in  Fignre  10),  which  is  most  closely  htted  aX  hm>  o.  as  essentially 


halfway  between  approximations  (Al)  and  (A2). 
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Fig.  1.  Regime  diagram  delineating  mesoscale  dynamical  responses  to  flow  U  across  the  elliptical 
obstacle  (1)  as  a  function  of  hm  and  (5,  following  Smith  (1989a, b).  Solid  circles  show  values 

whose  dynamics  we  model  in  this  study.  Dotted  gray  curves  show  critical  values  for  onset  of  wave 
breaking  {hm  >  hcb  =  F'<'cb)  and  flow  splitting  {hm  >  hcs  =  =  a). 
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(a) 


normalized  drag  D 


(b)  0  5  20 


Fig.  2.  Nondimensional  drag  for  obstacle  of  /3  =  1  and  hm  =  0.1:  (a)  vertical  profiles  of  at 
6,  12,  18  and  24  hours.  Gray  region  shows  sponge  layer.  Horizontal  dotted  line  shows  zt  =  6  km 
where  and  Dj  are  evaluated  and  plotted  in  (b)  as  time  series  together  with  surface  pressure 
drag  Dp  and  the  linear  t  >  oo  predictions  Dl  and  Dj. 
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Fig.  3.  Nondimensional  drag  as  a  function  of  time  for  the  (3  =  1  obstacle  for  normalized  obstacle 
heights  hm  ranging  from  (a)  0.1  to  (h)  1.5.  The  various  drag  curves  are  labeled  at  the  top  of  panel 
(a).  See  text  for  symbol  definitions. 
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Fig.  4.  Isentropes  as  a  function  of  x  and  z  (1  K  contour  interval)  along  the  central  obstacle 
axis  (y  =  0)  at  hourly  intervals  from  20-24  hours  for  hm  =  1-2  and  (3  =  1.  Color  shading  shows 
turbulent  kinetic  energy  (units  s“^:  see  color  bar).  Gray  line  depicts  Zi  =  6  km  where  gravity 
wave  momentum  fluxes  are  computed. 


42 


-20  0  20  40  60  80  -20 

X  (km) 

0 

20  40 

X  (km) 

60 

80 

-0.010  -0.008  -0.006  -0.004  -0.002  0.000 
kg  m'’  s'^ 

-0.001 

0.000 
kg  m"*  s’^ 

0.001 

-20  0  20  40  60  80  -20  0  20  40  60  80 

x(km)  x(km) 


-1.5  -1.0  -0.5  0.0  -0.3  -0.2  -0.1  0.0  0.1  0.2  0.3 

kg  m'’  kg  s'^ 


20  40 

X  (km) 


20  40 

X  (km) 


-2.5  -2.0  -1.5  -1.0  -0.5  0.0 

kg  m'’  s'^ 


-0.4  -0.2  0.0  0.2  0.4 

kg  m"*  s’^ 


20  40 

X  (km) 


20  40 

x(km) 


-3  -2  -1 

kg  m'’ 


-0.6  -0.4  -0.2  0.0  0.2  0.4 

kg  m’’'  s'^ 


Fig.  5.  Each  panel  pair  plots  the  FT-based  puw  and  pvw  fields  for  the  (5=1  obstacle  at  6  km 
altitude  after  24  hours,  for  hm  values  of  (a)  0.1,  (b)  1.0,  (c)  1.2,  and  (d)  1.5.  Flux  units  in  each 
accompanying  color  bar  are  kg  m~^  s“^. 
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Fig.  6.  Nondimensional  drag  as  a  function  of  time  for  the  /3  =  1  obstacle  for  normalized  obstacle 
heights  hm  ranging  from  (a)  2.0  to  (d)  5.0.  Curve  definitions  are  as  given  in  Figure  3a. 
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Fig.  7.  Isentropes  as  a  function  of  x  and  2  (1  K  contour  interval)  along  the  central  obstacle  axis 
(y  =  0)  at  hour  24  for  =  3  and  /?  =  1.  Color  shading  shows  turbulent  kinetic  energy  (units 
s“^:  see  color  bar).  Gray  line  depicts  Zi  =  Q  km  where  gravity  wave  momentum  fluxes  are 
computed. 
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Fig. 


As  in  Figure  5  but  for  hm  values  of  (a) 


2.0,  (b)  3.0, 


(c)  4.0,  and  (d)  5.0. 
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Fig.  9.  Nondimensional  drag  as  a  function  of  time  out  to  48  hours  for  the  /3  =  1  obstacle  for 
normalized  obstacle  heights  hm  of  (a)  1.2  and  (d)  3.0.  Curve  definitions  are  given  in  Figure  3a. 
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Fig.  10.  (a)  Nondimensional  and  (b)  absolute  mean  drag  {t  =20-24  hours)  and  range  (gray  error 
bars)  versus  hm  for  /3  =  1  WRF  simulations.  Dotted  gray  line  at  ~  1.1  demarks  responses  with 
and  without  TKE  at  smaller  and  larger  h^,  respectively.  The  various  drag  terms  are  labeled  top 
right:  see  text  for  symbol  definitions. 
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Fig.  11.  Nondimensional  drag  as  a  function  of  time  for  the  /3  =  3  obstacle  for  normalized  obstacle 
heights  hm  ranging  from  (a)  0.1  to  (h)  1.5.  The  various  drag  curves  are  labeled  at  the  top  of  panel 
(a).  See  text  for  symbol  definitions. 
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Fig.  12.  Isentropes  as  a  function  of  x  and  z  (1  K  contour  interval)  along  the  central  obstacle 
axis  {y  =  0)  at  hours  3,  5,  7,  9  and  11  for  hm  =  0.85  and  /3  =  3.  Color  shading  shows  areas  of 
turbulent  kinetic  energy  (units  s“^:  see  color  bar).  Gray  line  depicts  2*  =  6  km  where  gravity 
wave  momentum  fluxes  are  computed. 
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Fig.  13.  Nondimensional  drag  as  a  function  of  time  for  the  /3  =  3  obstacle  for  normalized  obstacle 
heights  hm  ranging  from  (a)  2.0  to  (d)  5.0.  Cnrve  definitions  are  as  given  in  Figure  11a. 
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Fig.  14.  (a)  Nondimensional  and  (b)  absolute  mean  drag  {t  =12-24  hours)  and  standard  deviations 
(gray  error  bars)  versus  hm  for  =  3  WRF  simulations.  Dotted  gray  line  at  hm  ~  0.7  demarks 
responses  with  and  without  TKE  at  smaller  and  larger  hm,  respectively.  The  various  drag  terms 
are  labeled  top  right:  see  text  for  symbol  definitions. 
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Fig.  15.  (a)  Nondimensional  mean  surface  pressure  drag  and  wave  drag  versus  hm  for  the  (5  = 
(solid)  and  (5  =  3  (dotted)  obstacles. 
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Fig.  16.  (Solid)  ratio  of  wave  drag  at  6  km  to  surface  pressure  drag  Dp  for  (a)  /3  =  1  and  (b) 
(5  =  3  obstacle.  Standard  deviations  in  (b)  were  derived  from  the  individual  standard  deviations  of 
Dw  and  Dp  in  Figure  14. 
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Fig.  17.  (a)  Normalized  and  (b)  absolute  mean  flow-blocking  drag  versus  hm  for  the  /3  =  1  (solid) 
and  (3  =  3  (dotted)  obstacles. 
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Fig.  18.  Schematic  diagram  of  the  mean  trends  (black  line)  and  range  of  variability  (thin  gray 
lines)  in  (a)  nondimensional  surface  pressure  drag  Dp  and  (b)  absolute  wave  drag  versus  hm, 
based  on  WRF  simulations  for  the  /3  =  1  and  (3  =  3  obstacles. 
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